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The detection of high-energy astrophysical neutrinos of extraterrestrial origin by the IceCube neutrino observatory in 
Antarctica has opened a unique window to the cosmos that may help to probe both the distant Universe and our cosmic 
backyard. The arrival directions of these high-energy events have been interpreted as uniformly distributed on the celestial 
sphere. Here, we revisit the topic of the putative isotropic angular distribution of these events applying Monte Carlo 
techniques to investigate a possible anisotropy. A modest evidence for anisotropy is found. An excess of events appears 
projected towards a section of the Local Void, where the density of galaxies with radial velocities below 3000 km s“^ is 
rather low, suggesting that this particular group of somewhat clustered sources are located either very close to the Milky 
Way or perhaps beyond 40 Mpc. The results of further analyses of the subsample of southern hemisphere events favour 
an origin at cosmological distances with the arrival directions of the events organized in a fractal-like structure. Although 
a small fraction of closer sources is possible, remote hierarchical structures appear to be the main source of these very 
energetic neutrinos. Some of the events may have their origin at the IBEX ribbon. 
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1 Introduction 

The IceCube neutrino observatory has announced the ob¬ 
servation of very high-energy events that have been inter¬ 
preted as the hrst robust evidence of neutrinos from ex¬ 
traterrestrial origin (Aartsen et al. 2013a,b,c, 2014a; Halzen 
2014) other than those detected from the supernova 1987A 
(see, e.g., Arnett et al. 1989) or the Sun (see, e.g., Bellerive 
2004; Bahcall 2005). A total of 37 events have been found 
in data collected from May 2010 to April 2013, an all-sky 
search over a period of 988 days, and their deposited en¬ 
ergies range from ^30 TeV to ~2 PeV. Assuming an as¬ 
trophysical origin for these detections, a number of possi¬ 
ble sources have been suggested (see, e.g., Murase et al. 
2013; Anchordoqui et al. 2014a; Joshi et al. 2014; Ahlers & 
Murase 2014). These include both Galactic and extragalac- 
tic sources (Halzen 2014) as well as others more exotic and 
rare (see, e.g., Bhattacharya et al. 2010). 

Among the proposed Galactic sources we hnd the cen¬ 
tre of the Galaxy (see, e.g., Razzaque 2013a), its halo (see, 
e.g., Taylor et al. 2014) with the Fermi Bubbles (see, e.g., 
Lunardini et al. 2014), and supernova remnants (see, e.g., 
Padovani & Resconi 2014). Suggested extragalactic sources 
range from active galactic nuclei (AGN; Winter 2013; KrauBI 
et al. 2014; Padovani & Resconi 2014), radio (Becker Tjus 
et al. 2014) or starburst galaxies (Anchordoqui et al. 2014b; 
Liu et al. 2014) to gamma-ray bursts (GRBs; see, e.g., Raz¬ 
zaque 2013b). Decay of a particle much heavier than PeV 
can be another possible source (Ema et al. 2014). Dark mat¬ 
ter annihilation has also been suggested as the origin of the 
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most energetic ones (Zavala 2014) although they can be ex¬ 
plained satisfactorily within the framework of the Standard 
Model of neutrino-nucleon interactions (Chen et al. 2014). 


Whatever the source, the angular distribution of events 
is considered consistent with an isotropic arrival direction of 
neutrinos (Aartsen et al. 2013c, 2014a) and the majority of 
events detected by IceCube are believed to be from cosmo¬ 
logical distances or, less likely, from the Galactic halo. In a 
recent update, Aartsen et al. (2014b) have found no evidence 
of neutrino emission from discrete or extended sources in 
four years of IceCube data. However, and although any dis¬ 
tribution of halo sources is expected to be isotropic, the 
spatial distribution of galaxies (star-forming or quiescent) 
located at cosmological distances is far from uniform, be¬ 
ing arranged in a hierarchical fashion (see, e.g., Pietronero 
1987). 


Here, we revisit the topic of the possible isotropic an¬ 
gular distribution of IceCube high-energy events using a 
Monte Carlo approach. This paper is organized as follows. 
The characteristics of the data used in our study are summa¬ 
rized in Sect. 2. Section 3 presents our Monte Carlo-based 
angular distribution analysis; an area of the sky is identi- 
hed as key to understand the most probable origin of the 
events studied here. The atmospheric background is further 
explored in Sect. 4. Our results are discussed within the con¬ 
text of the fractal structure of the Universe in Sect. 5. Con¬ 
clusions are summarized in Sect. 6. 


Copyright line will be provided by the publisher 




2 


de la Fuente Marcos & de la Fuente Marcos: On the angular distribution of IceCube events 



20 15 10 5 0 

at") 


Fig.l Distribution in equatorial coordinates (right ascen¬ 
sion and declination) of 36 events observed by IceCube. 
Data and errors as in Supplementary Table 1 of Aartsen et al. 
(2014a). Events with deposited energies in excess of 1 PeV 
(3) appear in blue (filled squares), the others in black (empty 
squares). Events containing muon tracks can be identified 
by their much smaller positional error bars. The Galac¬ 
tic plane/equator appears as a continuous black line. The 
Galactic disk is arbitrarily defined as the region confined 
between galactic latitude —5°and 5°(in green). The posi¬ 
tion of the Galactic centre is represented by a filled red cir¬ 
cle. The region enclosed between galactic latitude — 15°and 
15°appears in yellow, and the one confined by —30°and 
30°in cyan. See Table[T]for the number of events found to¬ 
wards these regions. The large gray empty circle signals the 
approximate location of the Local Void. Event 18 is the left¬ 
most one, event 10 is the rightmost one; both are well away 
from the bulk of the Galactic disk. Red empty circles signal 
probable background events (see the text for details). 


2 The data 


The data used in this study appear in Supplementary Table 
1 of Aartsen et al. (2014a). Thirty seven neutrino candidate 
events are included there, but event 32 was produced by 
a coincident pair of background muons from unrelated air 
showers and that made it impossible to reconstruct its arrival 
direction. This event will be excluded from our analysis as 
our initial working sample consists of 36 events (see Eig.[T]). 
Twenty-eight of them had shower-like topologies, while the 
remaining 8 had visible muon tracks. Three events have de¬ 
posited energies in excess of 1 PeV (in blue in Pig.[Tl)- Event 
14, with a deposited energy of nearly 1 PeV (the third high¬ 
est in the sample), is positionally coincidental (within the 
angular errors but almost 1?5 from Sagittarius A*) with the 
location of the supermassive black hole of the Milky Way. 
Twenty-seven events have negative declination and, there¬ 
fore, they arrived from the southern hemisphere. This sub¬ 
sample of 27 events is a time- and deposited energy-limited, 
and perhaps complete, sample. 




Fig. 2 Frequency distribution in equatorial coordinates 
(right ascension, top panel, and declination, bottom panel) 
of the 36 neutrino candidate events that appear in Supple¬ 
mentary Table 1 of Aartsen et al. (2014a). The number of 
bins is 2 where n = 36. An obvious excess of events 
in the right ascension range 8*^ to 20*^ is observed. 


It is widely accepted that the observed neutrinos come 
from two different sources, atmospheric and extraterrestrial. 
Atmospheric neutrinos are constantly being produced from 
all directions in the sky as they are the result of cosmic ray 
interactions in the upper atmosphere and subsequent decay 
of mesons; their deposited energies are expected to be un¬ 
der 100 TeV. The neutrino candidate events in Fig. [T] cor¬ 
respond to a deviation from the atmospheric background at 
the 5.1 (j level (Aartsen et al. 2014a). This means that out of 
37 events, about 11 are background and the remaining 26 
are extraterrestrial neutrinos. Aartsen et al. (2014a) already 
pointed out that the properties of events 3, 8, 18, 28 and 32 
are consistent with them being background detections. 

Figure |2] shows the frequency distribution in equatorial 
coordinates (a, 5) of the 36 neutrino candidate events dis¬ 
cussed above. There is an obvious anisotropy in declina¬ 
tion induced by the fact pointed out above that 27 events 
arrived from the southern celestial hemisphere. The source 
of this anisotropy is an observational bias, the IceCube neu¬ 
trino observatory is more sensitive to events arriving from 
that direction. There is also an obvious excess of events in 
the right ascension range 8^ to 20^. This range contains 27 
events, outnumbering by 18 the 9 events in the ranges (0, 
8)*^ plus (20, 24)^; this is a A.2a departure from an isotropic 
distribution, where a = ^Jn/2 is the standard deviation for 
binomial statistics and n is the number of events. It may 
be argued that the identified excess of neutrino candidate 
events towards the right ascension range 8** to 20** could be 
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an effect of the instrumental acceptance. However, the Earth 
rotation guarantees the uniform exposure of the detectors in 
right ascension (Aguilar et al. 2013). If we focus on the 27 
southern events, there are 18 of them within the right ascen¬ 
sion range 8^ to 20^ and this still translates into a somewhat 
significant 2.45(T departure from an isotropic distribution. 
Besides, this excess of neutrino candidate events towards 
the right ascension range 8** to 20** is consistent with the 
anisotropy in TeV cosmic rays identified by Schwadron et 
al. (2014). 

In their study, Schwadron et al. (2014) have found a 
broad relative excess of cosmic rays between a ~ 8^ and 
a ~ 18^. The details of the observed anisotropy depend 
on the cosmic ray energy range considered. The relative ex¬ 
cess in arrival direction is oriented towards the region of the 
sky that includes the downstream direction of the interstel¬ 
lar wind and the downfield direction of the Local Interstellar 
Magnetic Field (LIME; Schwadron et al. 2014). These TeV 
cosmic ray anisotropies are related to the Interstellar Bound¬ 
ary Explorer (IBEX) ribbon (McComas et al. 2009). Neu¬ 
trinos are not affected by the LIME, but high-energy cos¬ 
mic ray particles are and they can produce neutrinos when 
they interact with atomic nuclei creating unstable particles 
whose decay generates the neutrinos. As the parent high- 
energy cosmic ray particles can be deflected by magnetic 
fields, their original sources are almost impossible to trace. 

In this scenario, it is theoretically possible to produce ex¬ 
traterrestrial neutrinos at hundreds of AUs in the outskirts 
of the solar system; in fact, the structure observed by the 
IBEX spacecraft is a bright, narrow (~20°wide) ribbon of 
energetic neutral atom emission that likely hints at the pres¬ 
ence of an associated overdensity of matter. A putative pro¬ 
duction of neutrinos towards the edge of the solar system 
could further complicate the interpretation of the available 
data. 

Figure 0 shows the location of the IBEX ribbon struc¬ 
ture with respect to the Galactic equator and the IceCube 
events. The continuous yellow line follows the Local In¬ 
terstellar Medium (LISM) magnetic equator as derived by 
Heerikhuisen et al. (2010). A 20° wide band around the 
LISM magnetic equator appears in pink. The IBEX ribbon 
is found within that area. Some events appear to be aligned 
along the band. As the arrival distribution of high energy 
cosmic rays appears to be affected by perturbations induced 
in the local interstellar magnetic field by the heliosphere 
wake (see, e.g., Desiati & Lazarian 2012; Schwadron et al. 
2014) some of the detected neutrinos may have been gener¬ 
ated in that region. 

Returning to the connection between neutrino candidate 
events and discrete sources, Bai et al. (2014) have provided 
robust evidence on the possible Sagittarius A* origin of some| 
of the events observed by IceCube. In particular, event 25 
coincides spatially as well as temporally with some tran¬ 
sient activity observed in X-rays towards Sagittarius A* on 
2012 February 9. Another notable spatial coincidence is the 
galaxy 2dFGRS TGS070Z303 and event 18 (a muon-track 



R.A. (hours) 

Fig. 3 Distribution in equatorial coordinates (right ascen¬ 
sion and declination) of 36 events observed by IceCube as 
plotted in Fig.[T] The Galactic plane appears as a continuous 
black line. The position of the Galactic centre is represented 
by a filled red circle. The large gray empty circle signals the 
approximate location of the Local Void. The continuous yel¬ 
low line follows the LISM magnetic equator as derived by 
Heerikhuisen et al. (2010). A narrow ribbon (^20° wide) 
around the great circle defined by the LISM magnetic equa¬ 
tor appears in pink. This area matches approximately the 
location of the IBEX ribbon (McComas et al. 2009). 

event), separated by just over 59 arcseconds. Galaxy 2dF- 
GRS TGS070Z303 has a radial velocity of 37,866 km s"! 
(Colless et al. 2001). Also the galaxy 2dFGRS TGS280Z153| 
and event 10 (of shower-like topology) are separated by 
2(18. Galaxy 2dFGRS TGS280Z153 has a radial velocity 
of 53,795 km s“^ (Colless et al. 2001). An additional, per¬ 
haps significant, coincidence is the LINER AGN 2MASX 
J13550744-1309568 that is 5(16 from event 23, another muon-| 
track event at coordinates a = 13l‘91, 6 = —13?2. AGN 
2MASX J13550744-1309568 has a radial velocity of 24,924 km s" 
(Mauch & Sadler 2007). 

3 A Monte Carlo-based angular distribution 
analysis 

Let us consider a sample of n events observed projected 
on the surface of a sphere, the sky in our case. We want to 
study the possible non-uniform spatial distribution of these 
events. In order to do that, we generate random points on 
the surface of a sphere using an algorithm due to Marsaglia 
(1972). The positions of these points on the sky have right 
ascension in the range (O*', 24^) and declination in the range 
(-90°, 90°). We perform 10^ experiments, each one includ¬ 
ing n events, randomly distributed on the sky. In our his- 
tograms, the adopted number of bins is 2 rip , where rip 
is the number of pairs. We compute the angular distance, 
d, between two given points (or events) using the equation 
cosd = cos(ai — a 2 ) cosfji cos( 52 + sini5i sini52, where 
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(«!, (5i) and ( 0 : 2 , S 2 ) are the equatorial coordinates of both 
points. Then count the number of pairs per bin and divide 
by the number of experiments performed. 

The black continuous curves in Fig.|4]show the expected 
distribution of separations between events on the sky if they 
are uniformly distributed. The most probable separation is 
tt/2 radians or 90°. The distribution associated with the 
data in Supplementary Table 1 of Aartsen et al. (2014a) is 
far from regular (top panel) but only one bin deviates ^2cr 
(2.25cr) from a uniform distribution. Here, we use the ap¬ 
proximation given by Gehrels (1986) when rip < 21 to 
compute the error bars as Poisson statistics applies: a ^ 
1 + 1^0.75 + Tip (a = otherwise), where rip is the 
number of pairs. In principle, all the events are uncorrelated 
and distributed according to Poisson statistics. The addition 
of two independent Poissonian variables (considering ter¬ 
restrial and extraterrestrial events, see above) is still Poisso¬ 
nian. 

Figure IH top panel, does not take into account the fact 
that a fraction of the events (11 out of 37) likely have ter¬ 
restrial origin. Only the extraterrestrial neutrinos are of in¬ 
terest here, i.e. the sample must be decontaminated. Back¬ 
ground rejection is implemented following the conclusions 
obtained by Aartsen et al. (2014a), and also using an algo¬ 
rithm in which we assume that they will have deposited en¬ 
ergies under 100 TeV and that (approximately and based on 
Fig. 2 in Aartsen et al. 2014a) 8 background events have 
energies in the range 30-53.3 TeV, 2 in the range 53.3- 
76.7 TeV, and 1 in the range 76.7-100 TeV. In addition, we 
rank events based on how high is their probability of be¬ 
ing part of an isotropic distribution. At the energies of inter¬ 
est here, the angular distribution of atmospheric neutrinos is 
expected to be very nearly isotropic. We compute the sepa¬ 
ration distributions for all the possible combinations (more 
than 600 million) of 11 events and compare with that of 
an isotropic distribution of the same size. Those combina¬ 
tions of 11 events with root-mean-square deviation from the 
isotropic behaviour < 0.0047 are selected and the number 
of times each event appears is counted. Statistically speak¬ 
ing, those with the highest number of counts are more likely 
to be part of an isotropically distributed subsample. Using 
the conclusions from Aartsen et al. (2014a) and the events 
ranked as described, we conclude that the most probable 
background events are: 1, 3, 8, 9, 18, 27, 28, 29, 31, 32 and 
37. These events appear as red empty circles in Fig.[T](for a 
different approach to the issue of background rejection see 
Padovani & Resconi 2014). 

If we focus on the 26 events (or 325 pairs) which are 
most probably of extraterrestrial origin, we obtain the bot¬ 
tom panel in Fig. |4] For this decontaminated sample, the 
difference between the theoretical (isotropic) value and the 
actual one is always < 2.5 ct. However, the range 0 to 7r/2 ra¬ 
dians contains 195 pairs, outnumbering by 65 the 130 pairs 
in the range 7r/2 to tt radians; this is a 3.6 (t departure from 
an isotropic distribution, where a = ^Jn^j2 is now the stan¬ 
dard deviation for binomial statistics. Therefore, the irregu- 




Fig. 4 Separation distribution for the events in Supple¬ 
mentary Table 1 of Aartsen et al. (2014a), red line/points 
with error bars, the number of pairs is rip. The theoreti¬ 
cal separation distribution for an uniformly spread sample 
across the sky of n = 36 events (or 630 pairs, top panel) and 
n = 26 events (or 325 pairs, bottom panel) is displayed as 
a continuous black line. The top panel displays the entire 
sample of 36 events and the bottom panel shows the results 
excluding events of probable terrestrial origin, a subsample 
of 26. The blue lines correspond to the average separation 
distribution when the positional errors associated with the 
events are taken into account (see the text for details). The 

1(7 error bars have been computed as described in the text. 

1 /3 

The number of bins is 2 rip' . 

lar distribution on the celestial sphere is different from that 
of a sample of points uniformly scattered, i.e. an anisotropy 
exists. An excess of events appears to be located towards 
the equatorial coordinates (181*7, —23°); this has already 
been pointed out by Aartsen et al. (2013c). This region is 
relatively far from the Galactic disk, the densest part of the 
Milky Way Galaxy, and towards the Local Void (Tully & 
Fisher 1987). This area of the sky has a conhrmed low local 
density of galaxies with radial velocities below 3,000 km 
s“^ (see, e.g.. Fig. 2 in Nasonova cfe Karachentsev 2011). 
This region of low density extends up to distances of 40 to 
60 Mpc centred at least 23 Mpc from the Milky Way (Tully 
et al. 2008) and it may be even larger, perhaps connecting 
with the Microscopium/Sagittarius void (Kraan-Korteweg 
et al. 2008). The geometric centre of the Local Void is found 
towards the equatorial coordinates (19^*0, h- 3°), starting at 
the edge of the Local Group (Nasonova & Karachentsev 
2011). The large gray empty circle in Fig. [T] signals the 
approximate location and projected extension of the Local 
Void. Taking into account the error bars, it encloses nearly 
25% of all the IceCube events. 

On the other hand, the Galactic disk is the brightest 
source of gamma rays in the sky, can the excess of events 
pointed out above have their origin in the disk? To test this 
hypothesis we count the number of IceCube events within 
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Table 1 Number of events detected by IceCube as a func¬ 
tion of the galactic latitude, b, or distance from the Galactic 
plane in steps of 5°and the expected number of events de¬ 
rived from an isotropic model (see the text) with its asso¬ 
ciated error. The difference between the observational data 
and the isotropic model in units of cr is also given. Results 
for the decontaminated sample appear in boldface. See Fig. 
[T]for a reference on the various | A6| regions. 


lAb| n 

isotropic (dicr) 

IceCube data 

Difference (in units of a) 

5 

3.14±2.97 

3 

0.05 


2 . 27 ± 2.74 

2 

0.10 

10 

6.25±3.65 

10 

1.03 


4 . 51 ± 3.29 

8 

1.06 

15 

9.32±4.17 

15 

1.36 


6 . 73 ± 3.73 

11 

1.14 

20 

12.31±4.61 

16 

0.80 


8 . 89 ± 4.10 

12 

0.75 

25 

15.21±4.99 

18 

0.56 


10 . 98 ± 4.42 

14 

0.68 

30 

18.00±5.33 

20 

0.38 


13 . 00 ± 4.71 

16 

0.64 


bands whose boundaries are limited in galactic latitude, b 
(see Table [U. For example, points within the Galactic disk 
have |A&| = 5°and are confined between galactic latitude 
—5°and 5°(green area in Fig.[T]). We repeat the process us¬ 
ing the same uniformly distributed model studied before. 
Differences between the observational data and the isotropic 
model are not statistically significant. It appears unlikely 
that the dominant source of these high-energy neutrinos is in 
or near the Galactic disk. This is somewhat consistent with 
results obtained by Aartsen et al. (2014a); in their work, it 
is found that the correlation with the Galactic plane is not 
significant, with a probability of 2.8% for a width of the 
disk of ±7?5. Restricting the analysis to the range of right 
ascension 16*' to 22*“ does not show any significant changes 
and that also makes the Galactic bulge an unlikely source. 
Results are equivalent for the full sample and the decontam¬ 
inated one (see Table[il in boldface). 

Our previous analysis has not taken into consideration 
that the positional errors associated with the events, in par¬ 
ticular to those of shower topology, are very large, much 
larger than those typical of optical astrometry. What is the 
impact of these errors on our previous analysis? In order 
to evaluate the effect of large angular errors on the separa¬ 
tion distribution, we generate sets of synthetic events Monte 
Carlo-style (Metropolis & Ulam 1949) using the equato¬ 
rial coordinates and median angular errors in Supplemen¬ 
tary Table 1 of Aartsen et al. (2014a). For a given event of 
coordinates (a, (5) and median angular error Ae, we gen¬ 
erate a synthetic event of coordinates (as,<5s) where the 
new coordinates are such that the angular distance between 
the synthetic coordinates and the original ones is given by 
d* = |GRND| Ae, with GRND being a random number 
in the interval (-1, 1) following a Gaussian distribution ob¬ 
tained using the Box-Muller method (Press et al. 2007) with 


zero mean and unit variance. If the North Celestial Pole is 
within the error circle, declination is restricted to the range 
6 — Ae to +Tr/2, when the South Celestial Pole is included 
the declination range is —7r/2 to (5 -f Ae; in both cases, the 
right ascension can take any value in the range 0*“ to 24*“ as 
parts of all meridians are within the error circle as well. This 
approach is not singular at the poles. When none of the poles 
is included within the error circle, the range in declination 
is i5 — Ae to 5 -|- Ae and the equivalent in right ascension 
is a — Ae to a -f Ae. The average separation distribution 
for 10*^ sets of synthetic events based on the original ones 
and their positional errors is displayed as a blue line in Fig. 
|4](top panel). Taking into account the positional errors sig¬ 
nificantly reduces the difference between the observational 
angular separation distribution and that of an isotropic case. 
However, for the clean sample (bottom panel) there is still a 
3.1fT departure from an isotropic distribution when consid¬ 
ering separations smaller than 7r/2. The approach outlined 
here has also been used to strengthen the results of our de¬ 
contamination algorithm (see above) and our list of most 
probable background events has been obtained after hun¬ 
dreds of trials that use the median angular errors. 

4 The atmospheric background 

The background produced by atmospheric neutrinos is of 
diffuse nature. These atmospheric neutrinos result from cos¬ 
mic rays interacting with the atmosphere. So far, our analy¬ 
sis has assumed that this background follows a mathemati¬ 
cally isotropic model as the theory predicts. However, in the 
calculations carried out by the IceCube collaboration, the 
significance of their results is evaluated by comparing with 
a distribution obtained by performing yet another analysis 
over a sample of only-background data sets that are obtained 
by scrambling real data in right ascension. In other words, 
right ascension is randomized while keeping the original 
values in declination. This is done to take into account that 
the detector has anisotropic acceptance, i.e. its instrumen¬ 
tal response function is unknown. This approach assumes 
that the unknown instrumental response function does not 
change over time, i.e. it remains the same over the entire 
search period. Figure |5] shows that under this approach the 
separation distribution of the 27 southern events and that of 
an average of 10*^ only-background data sets are compati¬ 
ble, i.e. the angular distribution of the events is isotropic. 
However, deviations close to 1.5cr are observed. 

5 Discussion 

After being emitted, neutrinos move unaffected by gas, dust 
or magnetic fields. The somewhat significant excess observed! 
towards the equatorial coordinates (18(*7, —23°) is an im¬ 
portant clue to the origin of the observed events as their as¬ 
sociated neutrinos have travelled very little distance (if their 
sources are Galactic) or across mostly empty space in their 
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Fig. 5 Separation distribution (red line/points with error 
bars, top panel) for the 27 events located towards the south¬ 
ern hemisphere (negative declination) in Supplementary Ta¬ 
ble 1 of Aartsen et al. (2014a), the number of pairs is rip. 
The analogous separation distribution for an uniformly dis¬ 
tributed sample across a semisphere of n = 27 events or 531 
pairs is displayed as a continuous black line. The continu¬ 
ous green line shows the average of 10^ experiments where 
the right ascension is randomized as described in the text. 

final leg, the Local Void (if they are of extragalactic ori¬ 
gin). This is probably the key fact to understand the place of 
origin of the majority of events detected so far by the Ice- 
Cube neutrino observatory, extragalactic versus local. Our 
analysis shows no statistically significant excess of events 
associated with the Galactic disk (see above). If sources in 
the Galactic halo are expected to follow a uniform spread, 
then our findings are somewhat incompatible with such ori¬ 
gin. As for sources within the Local Group, multiple events 
are observed towards an area of the sky where the density of 
galaxies with radial velocities below 3,000 km s“^ is quite 
low. The alternative explanation is straightforward: sources 
located at cosmological distances, well beyond the Local 
Group. 

The spatial distribution of galaxies situated at cosmo¬ 
logical distances has been found to be hierarchical (see, e.g., 
Pietronero 1987) and this is also consistent with numerical 
expectations since early studies (e.g. Aarseth et al. 1979). 
In an attempt to further explore this interesting possibility, 
a dominant extragalactic origin for these events, we restrict 
our analysis to the sample of 27 events located towards the 
southern hemisphere in Supplementary Table 1 of Aartsen 
et al. (2014a). The separation distribution of this subsample 
(on the surface of a semisphere) is quite different from that 
in Fig.m with a reduced number of pairs separated by more 
than 7r/2 radians (the average angular distance for points 
uniformly spread on a sphere) with respect to the full sam¬ 
ple (see Fig.|6]l. 

In order to check if the angular distribution of events 
on the semisphere is isotropic, we proceed as we did with 
the full sample and perform 10^ experiments, each one in¬ 
cluding 27 events randomly distributed on the southern half 
of the celestial sphere, counting the number of pairs within 
a certain range of separations. The results are plotted as a 
black continuous line in Fig. |6] Differences between the 
uniform case and the subsample are consistent across the 
entire separation range and they are marginally significant 
for some bins (the largest difference is l.Scr). Our analysis 


shows that, for events observed towards the southern hemi¬ 
sphere, the angular distribution is not isotropic. There is a 
3.7(7 departure from an isotropic distribution when consid¬ 
ering separations smaller than 7r/2; as pointed out above, the 
deviation for the entire celestial sphere is 3.6(7. If we replace 
the original set of 27 data by sets of synthetic events ob¬ 
tained from the original ones and their Gaussian positional 
errors (see above) and average the results of 10^ experi¬ 
ments, we obtain the blue curve in Fig. |6] (top panel) that 
deviates even more clearly from the isotropic arrival sce¬ 
nario (black curve) for separations larger than the average 
(the largest difference is now close to 3a near the high end 
of the distribution). We observe 98 pairs with angular sep¬ 
aration in the range 7r/2 to tt instead of the expected 149 in 
the isotropic case, a 5.5cr deviation. This reduced number 
of pairs characterized by large separations strongly favours 
the existence of some type of weak but noticeable cluster¬ 
ing for these neutrino sources. Absence of strong anisotropy 
is usually interpreted as an indication of cosmological ori¬ 
gin. Using the decontaminated sample (see Fig. |6l bottom 
panel), we observe a 3.9 ct departure (149 observed pairs in¬ 
stead of 121 in the isotropic case for a total of 210 pairs) 
from an isotropic distribution when considering separations 
smaller than 7r/2. When taking into account the positional 
errors the statistical significance is reduced to 2.8a. 

If the sources of these neutrinos are located at cosmo¬ 
logical distances, probably within remote galaxy clusters, 
one may expect that they are organized following a hierar¬ 
chical fractal model. Star formation within the Milky Way 
follows a fractal pattern (e.g. de la Fuente Marcos & de la 
Fuente Marcos 2006a); similar patterns are observed in star- 
burst galaxies (e.g. de la Fuente Marcos & de la Fuente Mar¬ 
cos 2006b,c). Soneira & Peebles (1978) showed that the ap¬ 
pearance of the projected distribution of galaxies can be re¬ 
produced using a model based on the superposition of frac¬ 
tal clumps. In order to test this hipothesis we consider ria 
uniformly distributed points, (an, Sn), on the surface of the 
unit semisphere. Each point is the centre of an arc of radius 
i?/A. Within each arc we consider Ua additional uniformly 
distributed points, each one is centre for an associated arc of 
radius R/}?. Within each arc and on the surface of the unit 
semisphere, we generate uniformly distributed points. 
The approach is similar to the one described by Soneira & 
Peebles (1978). 

In our case, and in order to have 27 points (events) per 
experiment, we consider ria = 3, i? in the range ~0 to ~ 27r 
and A € (0, 2). After performing a large number of tests, 
each one consisting of 10® experiments, the best fit to the 
southern hemisphere data is found for R = tt and A = 0.9. 
The gray curve in Fig. |6] top panel, shows the distribution 
of separations for our fractal-like model. Although unable 
to reproduce the asymmetric behaviour at separations in the 
range 0.5-1.5 radians (top panel), its deviations from the 
actual data are not statistically significant and its overall 
agreement improves for the average separation distribution 
of synthetic events (blue curve in Fig. |6] top panel). These 
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Fig. 6 Separation distribution (red line/points with error 
bars, top panel) for the 27 events located towards the south¬ 
ern hemisphere (negative declination) in Supplementary Ta¬ 
ble 1 of Aartsen et al. (2014a), the number of pairs is Up. 
The analogous separation distribution for an uniformly dis¬ 
tributed sample across a semisphere of n = 27 events or 531 
pairs is displayed as a continuous black line in both pan¬ 
els. The average separation distribution for a set of synthetic 
events obtained from the original events and their Gaussian 
positional errors appears as a blue line, the average num¬ 
ber of synthetic events is ~25.8 and the associated number 
of pairs is ^319.9, more details in the text. The gray line 
shows the best fractal-like distribution (see the text for de¬ 
tails). The bottom panel plots the equivalent results for the 
decontaminated sample of 21 events. The 1 a error bars have 
been computed as described in the text. The number of bins 
is 2 Up , where rip is 531 (top) and 210 (bottom). 

deviations are always smaller than those from a uniform dis¬ 
tribution. In order to find out how significant the differences 
between the two models (fractal-like versus isotropic) are 
and because the sample is rather small, we use the likeli¬ 
hood formula due to Saha (1998). The average value of the 
W estimator for the fractal-like model is 6 times greater 
than the one for the isotropic model; therefore, it is statisti¬ 
cally more robust. The distribution of events over the south¬ 
ern hemisphere is almost certainly anisotropic. 

6 Conclusions 

In this paper, we have studied the angular distribution of 
IceCube high-energy events trying to answer the astrophys- 
ically relevant question of whether or not their spread is 
isotropic to conclude that the observed events are very prob¬ 
ably not evenly distributed on the sky. Differences between 
a uniform angular distribution and the observed one are in¬ 
deed marginally significant at the 2a level or better. Our 
main findings can be summarized as follows: 

- The angular distribution of IceCube high-energy events 
exhibits a modest but statistically significant anisotropy. 


- There is no statistically significant excess of events pro¬ 
jected towards the Galactic disk. 

- The observed distribution of events appears to be incom¬ 
patible with an origin in the Milky Way halo. 

- A relatively significant excess of events appears pro¬ 
jected towards the Local Void. This suggests either a 
Galactic origin or one well beyond the Local Group for 
these events. 

- The angular distribution of events towards the south¬ 
ern hemisphere appears to be anisotropic and it is sta¬ 
tistically more compatible with an origin in a distant, 
fractal-like organized structure than with the currently 
assumed isotropicity. 

We must stress that our results are based on small number 
statistics. However, they are compatible with well-establishecH 
results for GRBs, namely that short GRBs are distributed 
anisotropically (Balazs et al. 1998). Anisotropy is also ob¬ 
served for the highest energy cosmic rays (e.g. Abreu et al. 
2012; Schwadron et al. 2014). Because of the small size of 
the current data set compared to those of short GRBs or 
high-energy cosmic rays, the results of our study are neces¬ 
sarily limited in scope as only very large effects can be reli¬ 
ably identified. Although a small fraction of closer sources 
is possible (perhaps those putatively related to the IBEX rib¬ 
bon), remote hierarchical structures appear to be the main 
source of these energetic neutrinos. 
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